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Abstract 



We present a comprehensive investigation of light meson physics using maximally 
twisted mass fermions for Nj = 2 mass-degenerate quark flavours. By employing 
four values of the lattice spacing, spatial lattice extents ranging from 2.0 fm to 
2.5 fm and pseudo scalar masses in the range 280 < w-ps ^5 650 MeV we control the 
major systematic effects of our calculation. This enables us to confront our data 
with chiral perturbation theory and extract low energy constants of the effective 
chiral Lagrangian and derived quantities, such as the light quark mass, with high 
precision. 
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1 Introduction 



In recent years, the non-perturbative description of QCD on the lattice has 
made significant progress in tackling systematic effects present in the deter- 
mination of several important physical quantities (see e.g. [1,2] for recent re- 
views). Simulations containing the dynamics of the light-quark flavours in the 
sea, as well as those due to the strange quark and recently also to the charm, 
using pseudo scalar masses below 300 MeV, spatial lattice extents L > 2fm 
and lattice spacings smaller than 0.1 fm are presently being performed by sev- 
eral lattice groups. Such simulations will eventually allow for an extrapolation 
of the lattice data to the continuum limit and to the physical point while 
keeping also the finite volume effects under control. 

Adding a twisted mass term [3] to the Wilson-Dirac operator and tuning the 
theory to maximal twist [4, 5] has by now proved to be a practical and suc- 
cessful tool for performing (9(a)-improved lattice QCD simulations, see e.g. 
refs. [6-17]. A detailed study of the continuum-limit scaling in the quenched 
approximation [18-21] revealed not only that with an appropriate tuning pro- 
cedure to maximal twist lattice artefacts follow the expected [4] 0(a 2 ) scaling 
behaviour, but also that the remaining 0(a 2 ) effects are small, in agreement 
with the conclusions drawn in ref. [22]. The only exception seen so far is the 
neutral pseudo scalar mass [23] which shows significant 0(a 2 ) effects. We will 
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discuss this issue and its interpretation [24,25] further in sec. 3.4.1 below. 

In this paper, we shall present a study of the continuum limit scaling of lattice 
QCD with maximally twisted mass fermions for the case of dynamical fermions 
with Nf = 2 mass degenerate quarks demonstrating that also in this setup 
the 0(a 2 ) effects are small as has already been discussed in refs. [26-28]. 
We use data for the charged pseudo scalar meson mass mps and the decay 
constant /ps computed in a range of 280 MeV < mpg < 650 MeV. Our analysis 
is concentrated on two values of the inverse bare coupling 13,13 = 3.9 and 
P = 4.05 corresponding to values of the lattice spacing of a = 0.079 fm and a = 
0.063 fm, respectively. We complement our data set by adding an additional 
coarser lattice spacing of a = 0.100 fm corresponding to (3 = 3.8 and a finer 
lattice spacing of a = 0.051 fm, corresponding to f3 — 4.2. Our new results 
at (3 — 4.2, computed at two values of the pseudo scalar mass, confirm the 
smallness of lattice artefacts in the quantities considered here. 

The smallness of lattice artefacts, the range of masses covered, the usage of 
several lattice volumes and the precision with which mps and fps can be ob- 
tained using maximally twisted mass fermions allow us to confront the numer- 
ical data with chiral perturbation theory and to eventually extract accurate 
values for a number of low energy constants, in particular I 3 = 3.50(31) and 
£4 = 4.66(33), as will be discussed below. The main physical results we obtain 
from this analysis are the light quark mass mjff(/i = 2 GeV) = 3.54(26) MeV, 
the pseudo scalar decay constant in the chiral limit /o — 122(1) MeV, the scalar 
condensate [Y, m (pi = 2 GeV)] 1 / 3 = 270(7) MeV and / w //„ = 1.0755(94). The 
errors are statistical and systematical errors summed in quadrature. In case 
of asymmetric errors we use conservatively the maximum for the sum. The 
detailed error budget can be found in table 1 below. 

A special emphasis of this paper is the investigation of systematic errors on 
the results for the low energy constants from the above mentioned chiral fits. 
To estimate the systematic errors, we perform different kind of fits where 
we take lattice spacing, finite size and next-to-next-leading order (NNLO) 
corrections into account. Finally, we address the question of isospin violations 
by comparing the neutral to charged pseudo scalar masses as well as discussing 
the isospin splittings for other physical quantities. 

The paper is organised as follows: after describing the lattice set-up in the 
remainder of this section we shall present the main results of this paper in 
section 2. The following section 3 introduces the simulation set-up as well as 
the main ingredients entering in the continuum-limit scaling analysis of light 
meson observables. In section 4 we present the results of the combined chiral 
perturbation theory fits and section 5 gives a summary of our results. The 
details of the analyses and data tables are collected in the appendices. 
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1.1 Lattice Action 



The twisted mass action, the tuning procedure to maximal twist and the 
analysis techniques as employed by our European Twisted Mass Collaboration 
(ETMC) have been discussed extensively in refs. [6,9]. We therefore only 
briefly recapitulate the essential ingredients of our set-up. 

In the gauge sector we employ the tree-level Symanzik improved gauge action 
(tlSym) [29], viz. 

S 9 = 7;T,{ b ° E i 1 - B«^(0 + fe i E {! - ReTr(fO)^ 
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with the bare inverse gauge coupling (5 — 6/(?q, b± — — 1/12 and bo = 1 — 86i. 
The fermionic action for two flavours of twisted, mass degenerate quarks in 
the so called twisted basis [3,4] reads 

J2 fa) [D[U] +m + ipi^r 3 ] X (x)} , (1) 

X 

where m is the untwisted bare quark mass, /i q is the bare twisted quark mass, 
r 3 is the third Pauli matrix acting in flavour space and 



d[u] = \ [i, (v„ + v;) - av;v M 
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is the mass-less Wilson-Dirac operator. V M and V* are the forward and back- 
ward gauge covariant difference operators, respectively. Twisted mass fermions 
are said to be at maximal twist if the bare untwisted quark mass m is tuned 
to its critical value m crit , the situation we shall be interested in. For conve- 
nience we define the hopping parameter k — 1/(8 + 2am ). Note that we shall 
use the twisted basis throughout this paper. 

Maximally twisted mass fermions provide important advantages over Wilson's 
originally proposed formulation of lattice QCD: the spectrum of Q^Q with 
Q = j 5 (D[U] + m + ifiq'-fs) is bounded from below, which was the original 
reason to consider twisted mass fermions [3]. At maximal twist, the twisted 
quark mass \i q is related directly to the physical quark mass and renormalises 
multiplicatively only. Many mixings under renormalisation are expected to be 
simplified [4,5]. And - most importantly - as was first shown in ref. [4], physical 
observables are automatically 0(a) improved without the need to determine 
any operator-specific improvement coefficients. Another feature of maximally 
twisted mass fermions is that the pseudo scalar decay constant /pg does not 
need any renormalisation which allows for a very precise determination of this 
quantity. 
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The main drawback of maximally twisted mass fermions is that both parity 
and flavour symmetry are broken explicitly at non-zero values of the lattice 
spacing. However, it turns out that this is presumably only relevant for the 
mass of the neutral pseudo scalar meson (and kinematically related quantities). 
We shall discuss this issue in more detail later on. Note that in the following 
we use the notation mps = rrip S for the charged pion mass, while the neutral 
one will be denoted by mp S . 

Tuning to maximal twist is achieved with the general prescription [30, 31] 
to choose a parity odd operator O and determine am CTit such that O has 
vanishing expectation value at fixed physical situation for all lattice spacings, 
e.g. by fixing mps in physical units. One appropriate quantity to tune to zero 
is the PCAC quark mass, defined as 

Ex(«Wx,t)P"(0)) 

mpcAc= 2Ex(pa(X)t)pa(0)) , a = l,2, (2) 

where A ^ and P a are the axial vector current and the pseudo scalar density 
in the twisted basis, respectively, 

A U X ) = x(z)7m75 yX(z) , P a {x) = x(x)l5 yX(^) • 

Once rotated to the physical basis the numerator of the r.h.s of eq. (2) repre- 
sents the vacuum expectation value of a parity odd operator. Maximal twist is 
then achieved by demanding that the PCAC quark mass of eq. (2) vanishes. A 
discussion about the precise conditions we employ in our simulations is given 
in ref. [9] and summarised in section 3. 



2 Main Results 



In this section we present the main results of this paper. All the remaining 
sections are devoted to the discussion of the details leading to these results 
and to the estimation of systematic uncertainties. 

The ETM collaboration has gathered data for the pseudo scalar decay constant 
/ps and mass mpg with maximally twisted mass fermions at four different 
values of the lattice spacing ranging from 0.051 fm to 0.1 fm and for various 
values of the quark mass corresponding to values for the pseudo scalar mass 
from 280 MeV to 650 MeV. In addition to fp$ and mps we have determined the 
quark mass renormalisation constant = 1/Z P (at the three largest values 
of the lattice spacing only) and the hadronic length scale r /a. 

The Sommer parameter r [32] is used as a scaling variable in our investiga- 
tions, but not to set the physical scale. We have employed different ansatze for 
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(a) (b) 

Fig. 1. Scaling in finite, fixed volume for rj/pg at fixed values of r^mps (a) and for 
(rjmps) 2 at fixed values of r^/iR (b). In (b) we cannot include data at (3 = 4.2 due 
to the missing value of the renormalisation factor Zp. 

extrapolating r^/a to the chiral limit. While our data cannot discriminate be- 
tween these different ansatze, it turns out that all of our final physical results 
are independent of the exact way r /a is chirally extrapolated. This is due to 
the fact that the physical scale is eventually set by the pion decay constant 
f n . The values of r^/a we used for the continuum-limit scaling analysis were 
obtained by chirally extrapolating r /a linearly in (afi q ) 2 at every value of the 
lattice spacing separately, as explained in section 3. Later on, in the context of 
combined chiral fits, we also allow r /a to have an additional a\i q dependence, 
as discussed in section 4. 

Leading lattice artefacts in our data are expected to be of order a 2 , as discussed 
previously. This can be checked by extrapolating a physical quantity at fixed 
physical situation to the continuum limit. We show two such examples in 
figure 1. In the left panel we show / PS as a function of (a/rj) ) 2 at fixed value 
of r^mps- In order to match the values of r^mps at each value of r^/a and to 
fix the volume to ■ L = 5 we had to perform short inter- or extra-polations. 
The straight lines are linear fits in (a/rj) 2 to the corresponding data, with 
the data at the largest value of the lattice spacing not being included in the 
fit. It is clearly visible that the lattice artefacts appear to scale linearly in a 2 
and that their overall size is small. 

In the right panel we show the scaling of r^mps as a function of (a/r^) 2 at 
fixed values of the renormalised quark mass r^fi^, again at fixed, but finite 
volume. We conclude that also the charged pseudo scalar mass has only small 
lattice artefacts. This observation will become particularly important when 
we shall discuss later the mass difference between charged and neutral pseudo 
scalar mesons. Note that we do not yet have the value of Zp for the smallest 
value of the lattice spacing which is hence left out in fig. 1(b). 
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median 


statistical 


systematic 


m M [MeV] 


3.54 


(19) 


(+16 - 17) 


ft 

h 


3.50 


(9) 


(+9 - 30) 


h 


4.66 


(4) 


(+4 - 33) 


fo [MeV] 


121.5 


(0.1) 


(+1.1 - 0.1) 


B [MeV] 


2638 


(149) 


(±132) 


r r i 

r [fm] 


0.420 


(9) 


(+10 - 11) 


|S I 1 / 3 [MeV] 


270 


(5) 


(+3 - 4) 


U/ fo 


1.0755 


(6) 


(+8 - 94) 


/a(p = 3.90) 


5.32 


(5) 


/ i n\ 

(±0) 


r$/a(P = 4.05) 


6.66 


(6) 


(±0) 


a(/5 = 3.90) [fm] 


0.079 


(2) 


(±2) 


a(/3 = 4.05) [fm] 


0.063 


(1) 


(+1-2) 


Z P (/3 = 3.90) 


0.434 


(8) 


(+4-2) 


Z P (/3 = 4.05) 


0.452 


(9) 


(+3-9) 



Table 1 

Summary of fit results, determined from the weighted distribution as explained in 
the text. The first error is of statistical origin while the second, the asymmetric one, 
accounts for the systematic uncertainties. Bq, £ and m u d are renormalised in the 
MS scheme at the renormalisation scale [i = 2 GeV, as the values of Z-p are in the 
MS scheme at scale 2 GeV. Note that for the results listed here only data at j3 = 3.9 
and (3 = 4.05 have been used. The scale is set by = 130.7 MeV as done in ref. [6]. 
For a comparison to other recent lattice results we refer the reader to ref. [2]. 

The dependence of mps and fps on the renormalised quark mass and volume 
can be described by chiral perturbation theory (%PT) [33-35]. The residual 
lattice artefacts of order a 2 can also be included in the analysis. The corre- 
sponding formulae can be found in equation (8) below. We fit these formulae 
to our data in order to extract the parameters of the Nf = 2 chiral Lagrangian, 
i.e. the low energy constants and some derived quantities. Moreover, we can 
use these fits to calibrate our lattices by determining the value of the renor- 
malised quark mass Vq/z^ where the ratio mps//ps assumes its physical value 
(i.e. m^/fir) and set /ps = f n = 130.7 MeV there, [fj as done in ref. [6]. Hence, 
f n is used in this paper to set the scale. 

The results of these fits can be found in table 1. We give statistical and sys- 
tematic errors separately, the systematic one being asymmetrical. The results 
are obtained by performing O(80) fits, which differ in fit-range, finite size cor- 
rection formulae and in the order of xPT. The final result is obtained as the 
median of the corresponding weighted distribution over all fits. The statisti- 
cal error is determined using the bootstrap method with 1000 samples. The 
systematic uncertainty is estimated from the 68% confidence interval of the 
weighted distribution. For details we refer to the remaining sections of the 



2 Note that the recent update [36] for f n = 130.4(2) MeV does not affect our results 
within the errors we quote. 
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paper, in particular to appendix A. For a comparison to other lattice results 
we refer the reader to ref. [2]. 

The results quoted in table 1 are obtained from the two intermediate values 
of the lattice spacing only, since at the smallest available lattice spacing Zp 
has not yet been computed and the largest value of a is not considered, since 
there our condition for tuning to maximal twist might not have been fulfilled 
accurately enough, which may affect in particular the observables /ps and 
mps (at the smallest quark mass values) considered in this paper. However, 
including data at the smallest (by fitting = Zp 1 ) and the largest value of 
the lattice spacing gives results which are compatible within errors with the 
numbers shown in table 1, as can be seen from a direct comparison in table 7. 



3 Ensemble details, r , a 2 and finite size effects 

Before discussing the chiral perturbation theory fits to mps and /pg in section 4 
in detail we shall address in this section a number of issues which are needed as 
preparatory steps. These are the determination of the hadronic scale r /a, the 
tuning to maximal twist, the size of the 0(a 2 ) corrections (including isospin 
breaking effects) and finite size effects. 

3.1 Ensemble Details 

In table 2 we summarise the various Nf = 2 ensembles generated by the ETM 
collaboration. We remark here that all the configurations produced by ETMC 
are available on the international lattice data grid (ILDG) (see ref. [37] and 
references therein), and the configurations are made public after publication 
of this paper. We have simulations at four different values of the inverse gauge 
coupling (3 = 3.8, (3 = 3.9, (3 = 4.05 and (3 = 4.2. The corresponding values of 
the lattice spacing, obtained from chiral fits as detailed later, are a ~ 0.1 fm, 
a ps 0.079 fm, a ps 0.063 fm and a ps 0.051 fm, respectively. For each value of 
(3 we have several values of the bare twisted mass parameter a[i q , chosen such 
that the ensembles cover a range of pseudo scalar masses between 280 and 
650 MeV. 

The physical box length L of most of the simulations at f3 — 3.9, (3 = 4.05 are 
roughly equal and around L ps 2fm, while the volume at (3 — 3.8 is slightly 
larger. At f3 — 4.2, the lattice size ranges from L ps 1.7 fm to 2.5 fm. For 
all the values of /3, we demand the condition mpsL > 3. Note that we have 
carried out several simulations at different physical volumes for otherwise fixed 
parameters in order to study finite size effects (FSE). 
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Ensemble 


{L/af x T/a 


& 




K crit 


Tint(P) 


Tint(amps) 


T 


1 


24 3 x 48 


3.8 


0.0060 


0.164111 


190(44) 


8(2) 


1.0 


A 2 






0.0080 




172(80) 


10(2) 


1.0 


^3 






0.0110 




130(50) 


6(1) 


1.0 


A A 






0.0165 




40(12) 


6(1) 


1.0 


A 5 


20 3 x 48 




0.0060 




250(100) 


5(1) 


1.0 


Bi 


24 3 x 48 


3.9 


0.0040 


0.160856 


47(15) 


7(1) 


0.5 


B 2 






0.0064 




23(7) 


17(4) 


0.5 


B 3 






0.0085 




13(3) 


10(2) 


0.5 


B± 






0.0100 




15(4) 


7(2) 


0.5 


B 5 






0.0150 




30(8) 


20(6) 


0.5 


Be 


32 3 x 64 




0.0040 




37(11) 


2.8(3) 


0.5 


B 7 






0.0030 




51(19) 


7(1) 


1.0 


Ci 


32 3 x 64 


4.05 


0.0030 


0.157010 


18(4) 


7(1) 


0.5 


c 2 






0.0060 




10 2) 


9(2) 


0.5 


c 3 






0.0080 




13(3) 


7(1) 


0.5 


c 4 






0.0120 




5(1) 


4.8(6) 


0.5 


c 5 


24 3 x 48 




0.0060 




12(2) 


11(1) 


1.0 


c 6 


20 3 x 48 




0.0060 




10(2) 


7(1) 


1.0 


Dx 


48 3 x 96 


4.2 


0.0020 


0.154073 


13(2) 


< 8 


1.0 


D 2 


32 3 x 64 




0.0065 




6(1) 


< 8 


1.0 



Table 2 

Summary of ensembles generated by ETMC. We give the lattice volume L 3 x T and 
the values of the inverse coupling (3, the twisted mass parameter afj, g , the critical 
hopping parameter K cr it as determined at fJ, q , m m (i.e. the value of a\x q appearing in 
the same row as Kent) and the trajectory length r. The values of the lattice spacing 
that correspond to the four values of (3 are a ~ 0.1 fm ((3 = 3.8), a ~ 0.079 fm 
((3 = 3.9), a w 0.063 fm {(3 = 4.05) and a w 0.051 fm ((3 = 4.2). In addition we 
provide values for the integrated autocorrelation time of two typical quantities, the 
plaquette P and the pseudo scalar mass amps, in units of r = 0.5. We refer to 
ref. [9] for details on the determination of the autocorrelation time. 



The simulation algorithm used to generate these ensembles is a Hybrid Monte 
Carlo algorithm with multiple time scales and mass preconditioning. It is 
described in detail in ref. [38] and one implementation described in ref. [39] is 
freely available. Its performance for the current set-up is discussed in ref. [26]. 
In table 2 we provide the values of the actual trajectory length r used in 
the simulations. For each value of (5 and ajj, q we have produced around 5000 
equilibrated trajectories in units of r = 0.5. In all cases we allowed for at least 
1500 trajectories for equilibration (again in units of r = 0.5). Some of the 
ensembles have been produced with a trajectory length r = 1 following the 
suggestions of ref. [40]. 
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3.2 Hadronic scale parameter r 



In order to be able to compare results at different values of the lattice spac- 
ing it is convenient to use the hadronic scale r [32]. It is defined via the 
force between static quarks and can be measured to high accuracy in lattice 
QCD simulations. For details on our procedure to determine r /a we refer to 
ref. [9]. Note that in this work we use f n to set the scale and employ r only 
for a continuum limit scaling analysis. In other words is used only as an 
intermediate reference quantity, which is finally eliminated in favour of /„-. 

In fig. 2 we show r /a as a function of (afi q ) 2 for the examples of f3 — 3.9 and 
(3 = 4.05. The mass dependence appears to be rather weak and a linear fit in 
/iq describes the data well as has been also discussed in refs. [6,9]. However, 
also a linear dependence on a/i q cannot be excluded - due to possible effects 
of spontaneously broken chiral symmetry in tq. This dependence describes in 
the case of f3 — 3.90 the data even better than the quadratic ansatz. We shall 
therefore conservatively include both terms in our final analysis, i.e. the chiral 
fits discussed in section 4. 

For tuning to maximal twist as well as for the continuum-limit scaling analyses 
presented in this section - in particular in figures 1, 3 and 4 - we shall use 
the value of rj/a from a linear extrapolation in (afi q ) 2 to the chiral limit. A 
linear extrapolation in (afx q ) 2 yields r$/a = 4.47(6) at (3 — 3.8, r^/a = 5.25(2) 
at = 3.9, rlja = 6.61(2) at (3 = 4.05 and rf/a = 8.33(5) at (3 = 4.2. The 
statistical accuracy for r^/a is about 0.5%. It is important to notice that the 
ratio [ro/a(/3 = 3.9)]/[rg ja{J3 = 4.05)] is within the errors independent of the 
extrapolation procedure. This leads us to expect rather little influence of the 
r /a extrapolation strategy on physical quantities. 



3.3 Tuning to Maximal Twist 

Let us briefly discuss our strategy to tune to maximal twist. It corresponds 
to tune am to a critical value am crit such that the PCAC quark mass de- 
fined in eq. (2) vanishes. In particular, we determine the value of am CTlt at 
the lowest available value of an q>m i n <C oAqcd at each f3- value. Furthermore, 
we demand that the ratio |^AompcAc/ a A*g,mm| ^0.1 and also that its error 
A(ZAam PC Ac/a/ i g,min) % 0.1. These criteria and their justification for tuning 
to maximal twist have been discussed in ref. [9]. An additional condition is 
that the value of afi q>m i n corresponds to a fixed physical situation at all values 
of the lattice spacing. We have chosen to use mpg ~ 300MeV to define the 
value of a/ig jm i n . Considering (3 = 3.9, (3 = 4.05 and (3 = 4.2, it was possible to 
perform this tuning task with two or three tuning runs for each lattice spac- 
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0.0001 0.0002 0.0001 0.0002 



(ap,) 2 (afi q ) 2 

(a) (b) 

Fig. 2. ro/a as a function of (afi q ) 2 for (a) j3 = 3.9 and (b) j3 = 4.05. The lines 
represent a linear extrapolation in (an q ) 2 to the chiral limit. Note that we have 
always used the largest available volume for a given value of a/i q , see table 2. 

ing. Obeying these conditions leads to very small 0(a 2 ) effects in the physical 
observables considered in this paper, with the only exception of the neutral 
pseudo scalar mass, which will be discussed below. 

A special case are the simulations at (3 — 3.8. As can be seen in table 2 at the 
smallest value of a[i q = 0.006 we encounter very large autocorrelation times 
for the plaquette. Since the autocorrelation time for the PCAC quark mass is 
of the same order of magnitude as the one of the plaquette, at this value of a\i q 
it is questionable whether am-pcAc was reliably determined since the effective 
statistics is small. Hence, it is unclear whether a tuning to maximal twist 
fulfilling our conditions has been achieved at (3 = 3.8. A clarification of this 
point would have required substantially larger statistics. As a consequence, we 
will leave out the data at (3 = 3.8 in our main analysis. The (3 = 3.8 data will 
only be discussed in the context of systematic uncertainties. 

The measurements of mpcAc from all the ensembles are collected in Ap- 
pendix C. In figure 3 we show the properly renormalised ratio of the PCAC 
quark mass mpcAC over the twisted mass \x q (the renormalisation factor of Zp 
cancels in this ratio) against the renormalised twisted mass, //r, = fi q /Zp, in 
units of Vq, for the four values of f3. The renormalisation factors Zp, taken 
here at a scale of 2 GeV, were determined using the RI'-MOM scheme [41,42]. 
Note that for (3 = 4.2 the renormalisation factors are not available yet. How- 
ever, we have estimated their values by fitting Zp and by estimating the value 
of Za = 0.75 from the /3-dependence of its known values at smaller (3. Since 
these renormalisation constants enter as 0(1) multiplicative pre-factors, a 
small change in their values will not significantly change the data in figure 3. 
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■ (Z A m PC Ac/fJ-q) 


P = 3.80 

P = 3.90 i — ■ — i 

P = 4.05 • 

P = 4.20 i — ^ 






t. 





0.00 0.05 0.10 0.15 0.20 0.25 



Fig. 3. Renormalised ratio of the PCAC quark mass over the twisted mass against 
the renormalised twisted mass //r = Hq/Z-p at the four values of (3. The statistical 
uncertainties on Z-p and Z\ are not included. The data at j3 = 4.2 have been 
included by estimating the renormalisation constants as described in the text. At 
(3 = 3.8, the data for the lightest quark mass has not been included for the reasons 
explained in the text. The band indicates our condition for tuning to maximal twist, 
which is clearly achieved to a good precision. The arrow indicates the value of Tq^r 
where we tuned the PCAC mass to zero. 

Therefore our main conclusion about the quality of the tuning will be un- 
altered if the actual numbers for the renormalisation constants should come 
out slightly differently in the final analysis. 

As can be seen from fig. 3 our condition for tuning to maximal twist is satisfied 
with a good statistical precision. The fact that at the three largest values of (3 
the renormalised quark mass computed at the tuning twisted mass parameter 
a/x 9>min basically agree (this is illustrated by the arrow in fig. 3) shows in 
addition that the physical situation has indeed been fixed. At (3 — 3.8, due 
to the previously mentioned problem of the long autocorrelation time at the 
lowest quark mass, we include in fig. 3 only the ensembles for which a reliable 
estimate of the error on mpcAc was possible. At all (3- values for values of a\i q > 
o/%min we observe in ampcAc (small) deviations from zero. This 0(a) cut-off 
effect will modify only the 0(a 2 ) lattice artefacts of physical observables. 

3.4 0(a 2 ) effects 

In this section, we shall discuss the lattice spacing effects on the decay constant 
fps and the mass mps of the (charged) pseudo scalar meson. Starting with 
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/ps, we will use the data at (3 — 3.9, (3 = 4.05 and (3 = 4.2. In order to 
compare the data at these values of /3, we have fixed three reference values of 
the pseudo scalar mass to r^mps = 0.614, rjmps = 0.90 and r^mps = 1-10. 
The corresponding values for af PS were then obtained by small interpolations, 
which we have chosen to be linear in mp s . In addition, we also corrected for 
the very small differences in the physical volume for the (3 values used such 
that all data were scaled to the same physical volume • L = 5.0. This has 
been achieved by applying the relevant formulae from xPT as detailed below. 

As has already been discussed in section 2, in fig. 1(a) (cf. page 6) we show 
the pseudo scalar decay constant, brought to common reference points r^mps 
and volumes, as a function of (a/r^) 2 . In our earlier work, where we had only 
two values of the lattice spacing [26-28], we detected only very small lattice 
spacing artefacts when comparing (3 = 3.9 and (3 = 4.05. This observation is 
confirmed: in fig. 1(a) we now have added, for two values of afi q , a smaller 
lattice spacing of a = 0.051 fm and we can clearly observe that the scaling 
is as expected linear in a 2 with only a small slope such that no large cutoff 
effects are visible. Note that also for the nucleon mass the scaling violations 
are very small (see refs. [10,43]). 

We also show r^/ps for f3 = 3.8 in fig. 1(a) at least for those points where we 
are able to interpolate to the reference values of r^mpg. Although the linear 
continuum extrapolation is performed from the results at (3 — 3.9, (3 = 4.05 
and (3 = 4.2, as can be seen, the data at f3 — 3.8 are consistent with this linear 
behaviour. This indicates that the lattice artefacts are indeed small. 

In fig. 1(b) we show similarly the continuum extrapolation of (r^mps) 2 f° r 
three fixed values of the renormalised quark mass rg/iR = 0.045, 0.09, 0.13. In 
this figure we cannot include data at the smallest lattice spacing, due to the 
missing Zp-factor. Instead we include data from f3 = 3.80 in the continuum 
extrapolation to show that the scaling looks very flat. In conclusion also for 
(rg-mps) 2 residual scaling violations are small. 

3.4-1 Isospin Breaking Effects 

A peculiar 0(a 2 ) effect in the twisted mass formulation of lattice QCD is the 
breaking of isospin symmetry at any non-zero value of the lattice spacing. In- 
vestigations of these effects deserve therefore special attention. In the quenched 
approximation the difference (rj ) 2 ((mp S ) 2 — (m PS ) 2 ) between the neutral m PS 
and the charged rrip S pseudo scalar masses, although fully compatible with 
the expected 0(a 2 ) behaviour, turned out to be significant [23]. 

The determination of m PS requires the computation of disconnected diagrams 
rendering its determination difficult. Nevertheless, following the techniques 
described in ref. [9], we were able to compute m PS (and other neutral quan- 
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Ro 



m v /v 



fps 



3.90 
4.05 
3.90 
4.05 
3.90 
4.05 
3.90 
4.05 



0.004 
0.003 
0.004 
0.003 
0.004 
0.003 
0.004 
0.003 



0.04(06) 
-0.03(06) 

0.02(07) 
-0.05(09) 
-0.07(18) 
-0.31(29) 
0.022(29) 
-0.004(45) 



Table 3 



Comparison of some selected quantities for which an isospin splitting can occur 
for twisted mass fermions. Rq denotes the relative size of the splitting, see also 
ref. [24,25]. 

tities) with a reasonable precision. The results for m,p S and can be found 
in table CIO in the appendix. 

We show in fig. 4 the mass splitting between the charged and neutral pseudo 
scalar mesons, i.e. 



flavours the isospin breaking effects are large, although smaller than observed 
in the quenched approximation and of opposite sign, cf. ref. [23], meaning that 
in the case of dynamical quarks considered here the neutral pseudo scalar me- 
son turns out to be lighter than the charged one. This observation is consistent 
with the observed first order phase transition [44,45] and with the correspond- 
ing scenario of ref. [46]. 

One very important observation is that the large cutoff effect in the mass dif- 
ference between the neutral and charged pseudo scalar masses is dominated 
by discretisation effects in the neutral pseudo scalar mass. This is demon- 
strated in fig. 1(b) which shows the charged pseudo scalar mass as a function 
of (a/rj) 2 for fixed values of r^^ R . As can be seen, the lattice spacing effects 
are very small, as discussed above. 

Another observation is that, in the unitary sector, so far all other investigated 
differences between corresponding quantities affected by isospin violations are 
compatible with zero. In table 3 we have compiled the relative difference Ro = 
(O — 0')/0 for some selected observables and simulation points. Here O (O') 
denotes the charged (neutral) quantity in the case of mesons and A + (A ++ ) 
in the case of baryons. All of them are well compatible with zero. However, 
some quantities, like the vector meson decay constant fy, are rather noisy and 
hence the result for this quantity is not really conclusive. 

The ETM collaboration has been investigating the question of the large cut- 
off effects in the neutral pseudo-scalar mass [24,25,47] also theoretically. An 
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Fig. 4. The difference of the squared charged and neutral pseudo scalar masses 
as a function of a 2 in the Nf = 2 twisted mass formulation of lattice QCD at 
two different values of the charged pseudo scalar mass. A significant 0(a 2 ) lattice 
artefact is observed. The circles (triangles) correspond to a value of the charged 
pion mass of about 330 (430) MeV. The open circle is a larger physical volume. The 
lines are only to guide the eye and some points are slightly horizontally displaced 
for better visibility. Note that m PS was no ^ h e ld fixed for this plot, however, due 
to the large uncertainties on rrip S the picture should not significantly depend on 
this approximation. 



analysis a la Symanzik of the charged and the neutral pseudo scalar meson 
masses leads to the formulae 



Ks) 2 

(™ps) 2 



m 2 n + a 2 ^ + 0(a 2 m 2 n , a 4 
ml + 0(a 2 ml,a 4 ) 



(4) 



which show that the difference (mp S ) 2 — (^ps) 2 is given by the term propor- 
tional to £tt (for a discussion in the framework of twisted mass Wilson 
see ref. [48,49]). Here ( n = (ir°\CQ\7c ) and C§ is the dimension six term in the 
Symanzik effective Lagrangian. 



The main result of the analysis of ref. [24,25,47] is that is a large num- 
ber which in the vacuum saturation approximation can be estimated to be 
proportional to ICtt] 2 , where = (0|P 3 |7r°). The latter matrix element is 
numerically large: one finds \G n \ 2 / Aq CD around 20 — 25 [24,25,47]. This result 
provides a physical explanation for the large 0(a 2 ) effect in the neutral pseudo 
scalar mass. Moreover, since it can be shown that ( n appears only in the neu- 
tral pseudo scalar mass (and kinematically related quantities), one also finds 
a possible explanation of why all other splittings determined so far turn out to 
be small. We remark that in principle also a double insertion of the operator 
£5 can contribute to the 0(a 2 ) lattice artefacts. However, it turns out that 
these contributions are much smaller than those of £ 6 discussed here [25]. 
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3. 5 Finite Size Effects in /pg and mpg 



In the following, we shall focus on the treatment of finite size effects in /pg and 
m PS with the aim of assessing the applicability of an effective description to 
be used, in particular, to perform the very small volume corrections required 
to bring the simulation points to a common finite volume, as needed for the 
continuum-limit scaling analysis. Later on, in the context of the chiral fits 
described in section 4, we will apply the finite size corrections formulae to bring 
our continuum data to infinite volume and to check for systematic effects. 

Our lattice data have been obtained with values of mpgL > 3. It is generally 
believed that in such a situation the finite size effects are small and appear only 
at the per cent level. However, our simulation data are on a level of precision 
that such effects are certainly detectable and hence finite size effects can in 
general not be neglected. One possible way to control finite size effects is to 
use chiral perturbation theory (xPT) which provides an effective description 
for the dependence of physical quantities on the box size [50]. A completely 
independent analytical description of finite size effects (FSE) of particle masses 
is provided by the so called Luscher formula developed in ref. [51]. 

Being interested in mpg and /ps - where the applicability of %PT should be 
best justified - we then have two effective descriptions of FSE at hand: (i) 
xPT to next-to-leading order (NLO) by Gasser and Leutwyler [50], which 
we shall denote with GL for short, and (ii) the extended asymptotic Luscher 
formula of Colangelo, Diirr and Hafeli [52,53], which conveniently combines 
the Luscher formula approach with %PT and which we will refer to as CDH in 
the following. Note that in pure the finite size corrections to mps are also 
known to 2-loop order [54]. However, it turns out that the difference between 
the CDH formula and 2-loop %PT is negligible [54]. 

Below, we compare the numerical data obtained by ETMC to the two different 
analytical descriptions of the FSE in mpg and /pg, which is possible since we 
have precise values for these quantities for a number of ensembles such as C 2 , 
C 5 and Cq or B 1 and B 6 where all parameters but the lattice size are fixed. We 
remark that for all these ensembles mpgL > 3 holds, such that the formulae 
from refs. [50,52-54] should be applicable. This was also checked and discussed 
in ref. [26]. 

Starting with the GL ansatz, the finite size corrections for m PS and /pg read 



m PS (L) = m PS 1 + -&7i(A) 



m P g K m (L,m PS ,/o) , 
/ps Kf h (L,m PS ,f ), 



(5) 



/ps(£) = /ps 1-2&1 (A) 
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where 

f = mp s / (4vr/ ) 2 , A = m PS L , (6) 
gi is a known function [50] and the finite size corrections depend apart 
from L and mps only on the unknown leading order low energy constant fo 
representing the pseudo scalar decay constant in the chiral limit (note that 
our normalisation is such that the physical value of the pseudo scalar decay 
constant is f n = 130.7 MeV). 



The CDH formulae have a very similar form 



mps(L) = m PS 
fps(L) = fps 



n=l Z V nA 



= m PS K. 



CDH 

m 



70 m(n 



i + E7^{f+^+^? + o(rt} 

n=l V UA 



/ps ^ CDH 



/ 

(7) 



where the I^j can be written in terms of basic integrals as discussed in ref. [53], 
where also m(n) is tabulated. The authors of ref. [53] provide simplified for- 
mulae for the Imfi which we shall use in the following. Note that ij® is known 
while Iy is not. However, iffi is numerically so small that we shall drop it, 
i.e. set it to zero, in most of what follows, if^ = is always used. 



The drawback of the CDH formulae compared to GL is that additional low 
energy parameters are needed as an input[l], namely: Ai, A 2 , A 3 and A 4 . 
However, there are estimates for all those parameters given in ref. [55] and 
compiled in the form = log(A^/m^) in table 4 for convenience, which we 
shall rely on during the analysis presented in this section. In order to convert 
the estimates in table 4 into lattice units we conventionally used in this analysis 
= 0.45 fm. Note that the choice of scale here affects only the sub-leading 
contributions to the finite size corrections and has a negligible effect, as was 
checked also by direct inspection. We emphasise again that we are at this stage 
only interested in an effective description of FSE with the goal to understand 
whether the correction formulae are applicable in the regime of our simulation 
parameters. 

While the expansion in ref. [55] was performed in terms of the squared pion 
mass, in ref. [56] the same formulae have been obtained expanding in the 
quark mass directly. We shall refer to the latter as CDH m and include it in 
the following tests. 

We can now apply the different formulae to correct mps and fps for finite size 
effects for ensembles where we have different volumes available. The result 



3 If iffl is to be included there are also the parameters f needed. For estimates 
of their values see ref. [55] 
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i 




1 


-0.4 ±0.6 


2 


+4.3 ±0.1 


3 


±2.9 ± 2.4 


4 


±4.4 ±0.2 



Table 4 

Values for l- L = log(A?/m^) from table 2 of ref. [53] used here to analyse FSE. 





L/a 


amps(L) 


ampg (oo) 


CDH,6/ x 

am ps (oo) 


am^g H (oo) 


ampg Hm (oo) 


Bi 
Bq 


24 
32 


0.1362(7) 
0.1338(2) 


0.1354(7) 
0.1336(2) 


0.1348(7) 
0.1335(2) 


0.1350(7) 
0.1336(2) 


0.1346(7) 
0.1335(2) 


c 6 
c 5 
c 2 


20 
24 
32 


0.1520(15) 
0.1448(11) 
0.1432(06) 


0.1492(15) 
0.1436(11) 
0.1430(06) 


0.1448(15) 
0.1429(11) 
0.1429(06) 


0.1464(15) 
0.1432(11) 
0.1429(06) 


0.1430(15) 
0.1418(11) 
0.1427(06) 



Table 5 

Comparison of different finite size correction formulae of mps applied to ensembles 
where we have different volumes available, while all other parameters stay fixed. 
amps(L) is the measured value, ampg(oo) the infinite volume limit value from GL, 

arripg H,6 (oo) the infinite volume limit value from CDH, taking iffl into account and 
ampg H (oo) and ampg Hm (oo) the infinite volume limit value from CDH, without 

Im ■ In the last column the CDH formulae with a re-expansion in the quark mass 
is used. 





L/a 


afps(L) 


a/pf(oo) 


a/ P c s DH (oc) 


a/p C s DHm (-) 


Bi 
Bq 


24 
32 


0.0646(4) 
0.0663(2) 


0.0663(4) 
0.0666(2) 


0.0662(4) 
0.0665(2) 


0.0663(4) 
0.0665(2) 


c 6 
c 5 
c 2 


20 
24 
32 


0.0508(5) 
0.0558(5) 
0.0569(2) 


0.0546(5) 
0.0577(5) 
0.0573(2) 


0.0551(5) 
0.0574(5) 
0.0572(2) 


0.0558(5) 
0.0578(5) 
0.0572(2) 



Table 6 

same as table 5 but for fp$. 

can be found in tables 5 and 6. In the tables we show the corrected values at 
infinite volume from the different ways of treating finite size effects as indicated 
in the caption. To this end, we take the measured data at a given linear 
extent L and provide the corresponding correction to infinite volume. For 
larger volumes these extrapolations agree better and finally converge within 
the errors. By checking the stability of the values extrapolated to infinite 
volume from the different analytical formulae we can control the applicability 
of these corrections for a given lattice size. 

All finite size correction formulae provide an appropriate framework to de- 
scribe the observed finite volume effects. However, for ensembles B\ and Cq, 
corresponding to the smallest lattice extents at these two f3 values, GL seems to 
underestimate FSE in m PS , an observation which was also made in refs. [57,58]. 
For those cases the resummed Liischer formula provides a better description. 
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In the chiral fits presented in section 4, however, we observe significantly larger 
values of x 2 / d.o.f when GL is used, and hence we use CDH and CDH" 1 in these 
fits. 

In conclusion, these results make us confident that our simulations have even- 
tually reached a regime of pseudo scalar masses and lattice volumes where 
%PT formulae can be used to estimate FSE for m PS and / P s when ensembles 
with mpsL ^ 3 are used. In particular, this allows us to control the finite size 
effects for these mesonic quantities from our simulations. But it is clear that 
in particular the CDH formula is affected by large uncertainties, mainly stem- 
ming from the poorly known low energy constants, which are needed as input. 
Changing their values in the range suggested in ref. [53], however, changes 
the estimated finite size effects maximally at the order of about 20% (of the 
corrections themselves) [9]. 

Note that in the above discussion we have only used continuum formulae to 
describe the FSE of our lattice data. However, we believe that this can be 
justified by the smallness of the lattice artefacts discussed above. In addition, 
we restate that at this stage we are only interested in an effective but good 
description of the FSE since our goal was to bring the simulation data from 
different lattice spacings to a common finite volume as needed to perform the 
scaling analysis. In all the fits performed later on, the finite size corrections 
are applied to the continuum results. 



4 Combined Fits 



The main goal of this section is to confront the data for / PS and m P s to x?T 
and eventually determine low energy constants (LEC's) of the chiral effective 
Lagrangian. 

We shall proceed by presenting results from fits including a combined contin- 
uum, infinite- volume and chiral extrapolation of m PS and / P g for two /3-values, 
(3 = 3.9 and (3 = 4.05. Our analysis will extend the results given in refs. [26,27] 
by incorporating data for the renormalisation constant Z P and the Sommer 
parameter r^/a in the fit. Details on the computation of Z-p and r$/a can be 
found in refs. [9,42]. We shall include the correlation matrix of / P s and m P g 
into the fit (the correlations with r Q /a and Z-p are negligible) and also provide 
a detailed account of systematic uncertainties. 
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4-1 General Fit Formulae 



The formulae we are going to use to describe the chiral and continuum be- 
haviour of fps, mps in infinite volume read: 



(r£m P s) = (r£) x» 



1+^log 




-iNNLO 
■ m 



rofps = r$f 



l-^log (||)+ J D/ PS (a/rJ) 2 + T / NNLO 



(8) 



where T^j LO denote the continuum NNLO %PT terms 



T. 



NNLO 



E e 

2 S 



log 



jiNNLO _ _5^ 2 



log ir 



(9) 



with 



log 



A 2 

X\x 

M 



X„ 30 V " 
and we have defined 



Af 




A 2 


, A 2 


+ 32 log 




- 9 log ^ 


Xpt 






Xn 






A 2 


A 2 


+ 16 log 




+ 6 log 


Xn 




X/i 



A 2 
X^t 



(10) 



^Xm/(4^/o) 2 , XM = 2B /i R , » R = » q /Z? s (n = 2GeV). 



(11) 



We use a normalisation such that /o = V^-foj i- e - /tt — 130.7 MeV. For the 
finite size corrections we shall use the CDH and CDH m formulae as discussed 
in sub- sect ion 3.5. 

The mass and decay constant of the charged pion have been studied up to 
NLO [48, 49, 59] in the context of twisted mass chiral perturbation theory 
(tm^PT). The regime of quark masses and lattice spacings at which we have 
performed the simulations is such that /j,r > aAq CD . In the associated power 
counting, at maximal twist the NLO tnr\PT expressions for the charged pion 
mass and decay constant preserve their continuum form. The inclusion of the 
terms proportional to D mps j ps parametrising the lattice artifacts in eq. (8) 
represents an effective way of including sub-leading discretisation effects ap- 
pearing at NNLO. The data shown in fig. 1 (page 6) nicely support the /ir 
independence of -D/ PS ,m PS , which we assume in our fit ansatz (8). 



?o is the value of the Sommer parameter in the chiral limit which is determined 
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for every /3-value, within the global fit, according to the formula 



r /a) = (r*/a) 1 + C x r*fi R + C 2 (r^ R ) 5 



(12) 



with p independent, i.e. continuum parameters C\ and C 2 . In addition to the 
expected n R dependence, we also include the term linear in /i R in order to 
account for possible effects of spontaneously broken chiral symmetry. Possible 
lattice artifacts of 0(a 2 ) cannot be distinguished from the lattice artifacts 
included for r^mps and r^fps in eq. (8), and are therefore effectively described 
by the parameters D mps j PS , whereas lattice artifacts of order a 2 jiR and a 2 /z^ 
appear only at higher order. 

In order to work with a massless renormalisation scheme, the values of the 
pseudo scalar renormalisation constant Zp are determined in the chiral limit. 
They were determined using the RI'-MOM scheme at scale // = 1/a and 
the values can be found in table C.5. For the present fit all renormalisation 
constants Zp need to be run from scale fj,' — 1/a to a common scale, which we 
have conventionally chosen to be \x = 2 GeV. This running can be performed 
perturbatively, but the scale // is to be determined. This task can be performed 
for a given set of fit parameters by setting f n = 130.7 MeV where the ratio 
mps//ps assumes its physical value. Hence, we include for Zp 

Z]F(n = 2 GeV) = C(A*, AO Zf{fi' = 1/a) / R{fi = 2 GeV). (13) 

into the global fit. //) describes the running of Zp from scale // = 1/a 
to the scale /i — 2 GeV computed in perturbation theory to four loops in 
a s in ref. [60]. We use the value A^g = 0.2567 GeV [61]. The perturbative 
factor R(fi = 2 GeV) = 0.8178 matches the RI'-MOM to the MS scheme 
and is known to NNNLO. It should be noted that the scale \J enters only 
logarithmically into ( and the small changes in // therefore do make only a 
very small difference in the running of Zp. However, we include eq. (13) into 
the global fit in order to correctly account for the statistical uncertainty of 
Z«V = l/a). 

At NLO, i.e. setting T™ LO = 0, and neglecting finite size corrections for the 
moment, there are the following free parameters to be fitted to the data for 
a/ps, am PS , r /a and Zp 11 ' '(// = 1/a): 

r*/ , r*B , r *A 3 , 4 , C 1>2 , {r$/a} p , {Z^(2 GeV)}^, D mps , D fps , 

where we indicate with the notation {...}/3 that there is one parameter for each 
/3-value. When also the NNLO expressions for mps and fps are included there 
are four more parameters to be fitted to the data: 

^o A i,2, k M , k F . 
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As already mentioned, finite size effects are corrected for by using the asymp- 
totic formulae from CDH and CDH m , which is consistently included in the fit. 
However, both additionally depend on the low energy constants A 12 . 

It will turn out that we do not have sufficient data to determine all those 
parameters from the fit. For this reason we shall add priors for the NNLO low 
energy constants where needed. A description of the procedure to compute 
the x 2 can De found in appendix B. 

As mentioned above, in our analysis, we find that for a given value of afi q and 
f3 only the data for amps and a/ps are correlated and we use a correlated x 2 
to account for these. Estimates for the correlation between a/ps and amps can 
be found in the tables C.6-C.9 in appendix C. 

4-2 Analysis strategy 

By in- and excluding certain data points, we build different data-sets. For a 
given data-set we perform the following four different fits: 

(1) Fit A: NLO continuum xPT, T^JLO = 0j D mpsJps = 0, priors for r A lj2 

(2) Fit B: NLO continuum xPT, = 0, D mpsJps fitted, priors for r Ai )2 

(3) Fit C: NNLO continuum %PT, D mps j ps = 0, priors for roAi^ and ku,F 

(4) Fit D: NNLO continuum xPT, D mps j PS fitted, priors for r Ai i2 and Um,f 

In the above list of different fits, setting D mps j PS = corresponds to a con- 
stant continuum extrapolation of fp$ and mpg. All these fits are applied to 
every data-set summarised in table A.l, and are repeated with both CDH 
and CDH m , respectively. While we include or exclude m PS and fps data cor- 
responding to the data-set, for r /a we have always used only the largest 
volume available, since an effective description of finite volume effects is not 
available for this quantity. 

These different fits include the systematic uncertainties of lattice spacing arte- 
facts, finite size effects, higher order contributions in chiral perturbation theory 
and the extrapolation of r /a to the chiral limit. Our strategy is to build the 
distribution of all possible fits using the available data-sets, weighting with the 
corresponding confidence level. This procedure is analogous to what has been 
done in ref. [62]. In table 1 we give our best estimates for the fit parameters 
following this procedure. The best estimates and also quantities derived from 
those are obtained as the median of the weighted distribution from a total of 
76 different fits. All fits are repeated for 1000 bootstrap samples which are 
used to estimate the statistical errors, while the systematic uncertainties are 
estimated from the 68% confidence interval of the weighted distribution, as 
explained in appendix A. As an example we show in figure 5 the result for fit 
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Fig. 5. (a) (ro / a) / (tq / a) as a function of (r^/x^) 2 . (b) rj/ps, (c) (r^mps) 2 and 
(d) (r^mps) /(r^fin) as a function of r^p. The plots are for fit B on data-set 
1. Circles (triangles) represent data points from (3 = 3.90 {(3 = 4.05). The value of 
X 2 /dof obtained for this fit is 19/17. Note that in these figures we did not propagate 
the errors of Tq and Z-p. When these errors are included the statistical significance 
of a 2 -dependence of the data in the figures obviously decreases, see table 1. 



type B on data-set 1, which is also further discussed in the appendix. 

We exclude fits from the analysis which we consider as being not reasonable, 
i.e. where the NNLO fit gives a worse description at larger masses than the 
NLO fit. This is a sign for not having sufficient data for such a fit. In table 1 
there are a few quantities which show a strong asymmetry in the estimated 
systematics. This uncertainty is coming from the fact, that the estimate from 
NLO and NNLO fits for this quantities differs substantially. However, in the 
final result the NLO fits have more weight leading to the observed asymme- 
tries. 
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In order to obtain further confidence in the uncertainties we quote in table 1 we 
can use the data at /3 — 3.8 and (3 = 4.20. In table 7 we compare fits including 
these data with the results from table 1. It is clearly visible that the overall 
uncertainty we quote for all the quantities comfortably covers the differences 
among the three presented fit averages. This fact makes us confident that we 
indeed have a good estimate of the systematic uncertainties in our results. 
Note that Z P ((3 = 4.2) is a fit parameter for the fits including data at f3 — 4.2. 

We think that this strategy takes into account all systematic effects that may 
affect our data and which are accessible to us. This allows us to provide val- 
ues for the fit parameters and derived quantities with controlled statistical and 
systematic errors. In particular, since in our fit procedure, finite size correc- 
tions are applied in the continuum, effects of the neutral pseudo scalar mass 
on the finite size effects are absent. The details for our fits can be found in ap- 
pendix A. Note that the main uncertainty to our results is stemming from our 
estimates of systematic errors. A systematic effect we are not able to estimate 
is the missing strange and, possibly, charm quark in the sea. 

Eventually it is interesting to investigate the influence of the way r /a is 
chirally extrapolated on our final results. In table 1, the systematic effect 
on rj/a is negligible because of the little variation in the treatment of r /a 
between the different fits. Looking at table 7 one might be surprised that the 
values for C\ and C2 are both compatible with zero, even if in figure 5(a) a clear 
mass dependence is visible. The reason for this is that with both parameters 
the fit is overdetermined and the boostrap samples for C\ and C2 are strongly 
anticorrelated (—0.99). However, since we are in this analysis not primarily 
interested in the estimates for C± and C2, but more in the correct estimate of 
systematics for other quantities, we keep both parameters in the fits. 

In order to address this systematic effect, we have repeated all the analysis 
using 



i.e. setting C\ = 0, as compared to eq. (12). In table 8 we compare the results 
obtained with eq. (14) with the ones presented in table 1. The comparison 
leads to the conclusion that only the values of rj/a are significantly changed, 
being also statistically more precise when eq. (14) is used. Thus, the system- 
atic uncertainty stemming from the r^/a extrapolation to the chiral limit is 
negligible for physical quantities. 

In the tables 1, 7 and 8 we also quote a number of quantities that are derived 
from the basic fit parameters. The values for the low energy constants £j, 



(r /a) = (rgf/a) 1 + C 2 {r^ R f 



(14) 




(15) 
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/3 = 3.8, 3.9, 4.05 


(5 = 3.9,4.05 


(3 = 3.9, 4.05, 4.2 


prior 


m ud [MeVl 


3.84(18) 


3.54(26) 


3.58(26) 




h 


3.32(21) 


3.50(31) 


3.49(9) 






4.69(17) 


4.66(33) 


4.63(4) 




h 


-0.41(58) 


-0.33(62) 


-0.59(58) 


-0.4(6) 


h 


4.31(11) 


4.32(11) 


4.32(10) 


4.3(1) 


fo [MeV] 


121.7(3) 


121.5(1.1) 


121.58(8) 




B [MeV] 

w L J 


2437(120) 


2638(200) 


2619(190) 




r [fm] 


0.446(9) 


0.420(14) 


0.429(8) 




Ci 


-0.09(18) 


-0.21(17) 


-0.37(14) 




c 2 


-1.03(81) 


-0.52(77) 


0.09(60) 


_ 


/ r2) NLO r fm 2] 

\ / S L J 

lEl^s [MeV] 

II L J 


0.710(28 


0.715(77) 


0.710(9) 




262.2(4.0) 


269.9(6.5) 


268.4(6.6) 




U/fo 

r?/a(/9 = 3.80) 


1.0742(81) 


1.0755(94) 


1.0750(8) 




4.462 45) 








r?/a(/3 = 3.90) 


5.259(48) 


5.316(49) 


5.361(39) 




r?/a(/J = 4.05) 


6.637(57) 


6.661 62 


6.727(48) 




r*/a(/3 = 4.20) 






8.358 63 




a(8 = 3.80) ffiml 

V" / L J 


0.0998(19) 








a(P = 3.90) [fm] 


0.0847(15) 


0.0790(26) 


0.0801(14) 




a(/3 = 4.05 fm 


0.0672(12) 


0.0630(20) 


0.0638(10) 




a(/3 = 4.20) [fm] 






0.05142(83) 




Z P (0 = 3.80) 


0.431(10) 









Z P (/3 = 3.90) 


0.4390(74) 


0.4335(84) 


0.436(9) 




Z P (/3 = 4.05) 


0.455(10) 


0.452(13) 


0.450(11) 




Zp(P = 4.20) 






0.466(19) 


0.45(10) 




0(1) 


-0.7(1.4) 


0(2) 






0.56(45) 


1.68(68) 


1.5(6) 





Table 7 

We compare fit results of fits including data from (3 = 3.8 (first data column) and 
data from [3 = 4.20 (third data column) with the result where only data from 
(3 = 3.90 and (3 = 4.05 have been used (second data column). Note that in the 
case of including data from (3 = 4.20 the value of Zp^ s is determined from the fit. 
All averages are weighted averages with the confidence levels of the individual fits. 
The errors are statistical and systematical, added in quadrature. Bq, S and m u ^ 
are renormalised in the MS scheme at the renormalisation scale fi = 2 GeV, as the 
values of Zp are in the MS scheme at this scale [i. The scale is set by f n = 130.7 MeV 
as done in ref. [6]. 



are determined using the physical value of the charged pion mass m n , which 
is the usual convention. The scalar condensate can be determined from the fit 
parameters B Q and fo in the following way: 



B ofo 



(16) 
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table 1 



quadratic 



m u ,d [MeV] 


3.54(26) 


3.56(29) 


h 


3.50(31) 


3.52(48) 


h 


4.66(33) 


4.67(25) 


h 


-0.33(62) 


-0.50(68) 




4.32(11) 


4.33(11) 


fo [MeV] 


121.5(1.1) 


121.48(89) 


B [MeV] 


2638(200) 


2618(200) 


r [fm] 


0.420(14) 


0.416(14) 


Ci 


-0.21(17) 


— 


c 2 


-0.52(77) 


-1.45(14) 


(r 2 )NLO [ fm 2] 


0.715(77) 


0.720(60) 


IS] 1 / 3 [MeV] 


269.9(6.5) 


268.8(6.2) 


fn/fo 

r*/a{f5 = 3.90) 
r*/a((3 = 4.05) 


1.0755(94) 


1.0759(80) 


5.316(49) 


5.261(14) 


6.661(62) 


6.595(23) 


a((3 = 3.90) [fm] 


0.0790(26) 


0.0790(27) 


a((3 = 4.05) [fm] 


0.0630(20) 


0.0631(21) 


Z P (P = 3.90) 


0.4335(84) 


0.4341(90) 


Z P (/3 = 4.05) 


0.452(13) 


0.451(13) 




-0.7(1.4) 


0.1(1.7) 


Df PS 


1.68(68) 


2.10(87) 


Table 8 



We compare the main results of table 1, where the chiral extrapolation of tq was 
performed according to eq. (12), with results obtained with the purely quadratic 
extrapolation of r^/a to the chiral limit in eq. (14). For further details see table 1. 



and the scalar radius at next to leading order is given by 

(r 2 r° = TT^ih -13/12). (17) 



5 Summary and Outlook 



In this paper we have demonstrated that Wilson twisted mass fermions, when 
tuned to maximal twist, show the expected linear scaling in a 2 with an almost 
negligible coefficient for both the pseudo scalar decay constant / P g and mass 
mps (see fig. 1 on page 6). A very similar behaviour is seen for the nucleon 
mass [10,43]. It is important to remark that tuning to maximal twist has been 
achieved in our work by demanding that the renormalised ratio Zpjn-pcAc I A*™ m 
of the PCAC mass mpcAc over the minimal twisted mass /i™ m employed in 
our simulations satisfies | Zaitipqac/ ^ m I < 0.1. The same condition is used 
for the corresponding error of this ratio. 

The precise data for / PS and m PS we have obtained at several values of the 
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lattice spacing, quark masses and volumes allowed for a detailed analysis of 
systematic effects originating from having to perform a continuum, infinite 
volume and chiral limit. In particular, in this way it became possible to con- 
front our data to analytical predictions from chiral perturbation theory and 
to extract a number of low energy constants and derived quantities as listed 
in table 1, page 7. The main physical results we obtain from this analysis are 
the light quark mass m^(/i = 2 GeV) = 3.54(26) MeV, the pseudo scalar 
decay constant in the chiral limit /o = 122(1) MeV, the scalar condensate 
\E m (n = 2 GeV)] 1 / 3 = 270(7) MeV and U/f Q = 1.0755(94). 

By performing O(80) different fits, taking into account or leaving out terms 
describing e.g. lattice artefacts or NNLO chiral perturbation theory and in- 
cluding or excluding data sets, we could achieve a reliable estimate of sys- 
tematic errors. The final numbers for the fit parameters and their errors were 
then derived from the median, the bootstrap samples and the 68% confidence 
level of the (weighted) distribution of fit parameters originating from the dif- 
ferent fits. Some quantities such as Is = 3.50(31), U = 4.66(33) or the scalar 
condensate belong to the most precise determinations available today. For a 
comparison to other lattice results we refer the reader to ref. [2,63,64]. 

As a particular effect appearing in the twisted mass formulation of lattice 
QCD, we determined the size of isospin violation in various quantities. In this 
analysis we found, in accordance with theoretical expectations [24], that only 
the neutral pseudo scalar mass is significantly affected by the isospin breaking. 

We conclude that Wilson twisted mass fermions at maximal twist provide a 
lattice QCD formulation allowing for precise computations of quantities in the 
light meson sector of QCD. Clearly, the next step is now to include dynamical 
strange and charm degrees of freedom in our simulations and work in this 
direction is in progress [65,66]. 
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Appendices 

A Details of fits and discussion of systematics 

In this appendix we discuss the details of our data analysis. For a more com- 
prehensive discussion see ref. [36]. For every fit we compute besides the usual 
X 2 value also the associated confidence levels (CL) (sometimes called goodness 
of the fit) 

CL(g = X \n = dof) = 1 - P(q/2, n/2) , 

where P is the incomplete Gamma function which corresponds to the cu- 
mulative ^-distribution function for the given number of degrees of freedom 
F x 2(q,n). F x 2(q,n) itself represents (in the frequentists approach to statistics) 
the probability of finding a value of x 2 < Q- A confidence level of CL(g) = x 
indicates then the fraction of times the computed value of x 2 is larger than q, 
even if the fitted model is correct. 

Since we follow a Bayesian approach by including prior knowledge for some 
parameters into our fits, the interpretation of CL is more a Bayesian credibility 
level. Hence, we are not going to determine the most probably correct model 
on the basis of the values of the CL's, but we shall rather incorporate all 
different fits into the final result. The CL's of each fit are then used as weights, 
specifying the impact of a given fit on the final result. 

In more detail, for obtaining the final results we first sort out all physically un- 
reasonable fits, as will be described below. Then we generate the CL-weighted 
distribution of a given fit parameter or of a derived quantity 6 over all the 
retained fits. The expectation value of is estimated as the median of the 
weighted distribution. By performing this procedure for every bootstrap sam- 
ple we determine the statistical error on our estimate for 6. An estimate for 
the systematic uncertainty is provided by 68.4% confidence interval of the 
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weighted distribution. The estimate for the final error is obtained by adding 
the statistical and the systematic error in quadrature. 

These systematic errors cover effects from lattice artefacts, from different or- 
ders in xPT, finite size effects and fitting range of the quark masses. We present 
separately the results including data generated at /3 — 3.8 and f3 = 4.2, be- 
cause they are on a different level of accuracy due to insufficient precision in 
tuning to maximal twist at (3 — 3.8 and due to the missing estimate of Zp 1 
at (3 = 4.2. The results for these averages can be found in table 7, in- and 
excluding j3 — 3.8 and (3 = 4.2, compared to the result using only (3 — 3.9 and 
(3 = 4.05 data. The total residual error we quote for our final results covers 
the difference to the results including (3 = 3.8 and (3 = 4.2 data. 

As mentioned before, we vary the data-sets in order to probe the influence of 
the fit range and the finite volume effects on the fit. The data-sets we 
use are compiled in table A.l. We use both CDH and CDH m to correct for 
finite size effects to further check for the influence of those on our results. 

The priors for roAi^ and JtM,F are necessary to obtain stable fits. The priors 
for r Ai i2 are taken from table 4 and /cm,f = ± 10 is sufficient. 

NNLO Fits 

The NNLO %PT formulae provide in general a better x 2 /d.o.f. than the cor- 
responding NLO formulae. However, it seems that the larger masses with 
r oA t R ~ 0-14 are necessary to stabilise these fits. If these masses are left out 
from the fits then the NNLO fit curves strongly at large values of r^/iR and 
actually undershoots the data at larger masses. This can be seen in figure A.l, 
where we show in the left panel the data for r /ps as a function of rQ/x R and 
the best fit function for fit C. Note that the data at large values of rJ/i R are 
not included in the fit. In the right panel we show the LO, NLO and the com- 
plete NNLO function for comparison. Evidently, for larger masses, there is a 
significant difference between the NLO and NNLO results. It is in particular 
suspicious that the NNLO fit shows a stronger curvature than the NLO fit 
and that the NNLO fit deviates more from the data points at large masses 
- which are not included in the fit. On the contrary one would expect the 
NNLO formulae to describe the data points at larger masses better than the 
NLO formulae. On the other hand, when including the heavier masses, the 
NNLO fit is able to describe these data points. To improve the sensitivity of 
our lattice data to %PT at NNLO, additional data points would be needed. 

We remark that in ref. [56] the inclusion of NNLO was found to be 

necessary in order to obtain an adequate fit to the data on the pion radius. In 
ref. [56] it was argued that due to the p-meson dominance one should expect 
a slower rate of convergence of chiral perturbation theory. However, using the 
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Table A.l 

We list the various sets of ensembles used for the fits. We refer to table 2 for the 
specification of the ensembles used. 



data for a/ps and amps alone turns out to be insufficient to obtain the low 
energy constants at NNLO. 

We observe this behaviour for all NNLO fits for data sets without the ensem- 
bles at large values of Tq(x r > 0.14, therefore, we do not consider these fits 
as reasonable and exclude them in the weighted average, even if the CL for 
this fit would justify the inclusion. Note however that including them does not 
change the final results within the quoted uncertainties. 

As already discussed in section 4, one effect of the NNLO fits included in the 
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Fig. A.l. We show fit C for data set 1, i.e. including NNLO terms, but excluding 
lattice artifacts. This fit was not included in the weighted average. In (a) the best 
fit function for /ps as a function of the renormalised quark mass. In (b) we plot 
/ps//o i n LO, NLO and NNLO for comparison. 



final analysis is that the estimated systematic uncertainty in e.g. / , £3 and £4 
is asymmetric. If we analyse NNLO fits only, we would compared to table 1 
rather obtain f = 122.6(1.1), £ 3 = 3.22(66), £ 4 = 4.33(34). This uncertainty 
is included in our final error estimate of those quantities. 



Lattice Artifacts 

When performing fit B, which uses only NLO and includes lattice arti- 
facts, we find that the fitting parameter D mps parametrising the lattice arti- 
facts in mps is compatible with zero within errors. In /pg lattice artifacts are 
more apparent, even if the error on Df ps is still large. 

However, from the fits we clearly observe that the x 2 values are significantly 
smaller when lattice artifacts are included, and the weighted distribution is 
hence dominated by the fits of type B and D. As an example for such a fit we 
show in figure 5 plots for fit B on data-set 1. We show the data together with 
the fitted curves for {tq/cl)/{tq/ol) as a function of (ro/ipj 2 in (a), for rj/ps 
in (b), for (r^mps) 2 in (c) and for (ro?nps) 2 /( r o/ i R) m (°0, the last three as a 
function of rp/UR. Note that in figures 5(b), 5(c) and 5(d) we did not include 
the error of Tq and Zp for the data points. 

From figure 5 one can observe that the fit works very well, leading to % 2 /dof = 
19/17. In (a) one can see that the data for r Q /a is also fully compatible with 
a linear dependence on (rg/ipj 2 only, as was also discussed earlier and sum- 
marised in table 8. 
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Finite Size Corrections and Effects of Priors 

We have used the CDH and CDH m formulae for the fits and include both 
into the weighted distribution of all fits. It is still interesting to compare two 
fits, using CDH m on the one hand and the parameter free GL formulae on 
the other hand, while all the rest is identical. We observe that in general the 
% 2 -value using CDH m is about a factor of 2 smaller than the one using GL. 
Therefore, we think that having extra parameters (which are included in the 
fit with priors) is justified. 

Another test for uncertainties stemming from FS corrections is to perform 
fits where only the largest available volume per /i^-value is used. Including 
these ensembles (which correspond to data-sets 6,7 and 8) tend in general to 
produce improved values of x 2 , while the influence on other fit parameters is 
again small. Clearly, it would be desirable to have for all ensembles values of 
mpsL > 4, where all the FS correction formulae give very similar results and 
the uncertainty from finite volumes is hence much reduced. In this investiga- 
tion we include this significant uncertainty in our final error estimate. 

As discussed previously, we have to add priors for £ lj2 and /cm,f i n order to 
stabilise the fits. It is therefore interesting to discuss the influence of the choice 
for the priors on the fit results. First of all, the influence of the choice for Um,f 
and their errors has negligible effect on our fit results. 

In case of £i 5 2 the influence depends on whether or not NNLO terms are 
included in the chiral expansion. While doubling the error estimates of or 
changing their values within the error range has negligible effects in the case 
of a NLO chiral fit (fits A and B), many fit results change drastically in the 
case of NNLO chiral fits (fits C and D): we observe in the latter case that the 
best fit values e.g. for £3 and £4 differ substantially from the NLO fit result and 
also from the estimates given in table 4. At the same time the errors become 
much larger on those quantities, up to 50%. 

Moreover, the x 2 minimisation process becomes more difficult, the \ 2 functions 
shows many nearby local minima in parameter space. These facts confirm once 
more that /ps and mps data are not sensitive to NNLO terms in the chiral 
expansion. 

Chiral Extrapolation of Zp 1 ' 

The chiral extrapolation of Zp 1 ' is not performed in the fits, but the chiral 
value plus error estimates is input to the fit. In order to check whether this 
might lead the fits into a wrong direction, we have performed fits with doubled 
uncertainty on the values of Zp 1 . It turns out that this has negligible effect 
on our results. 



32 



B x 2 Implementation 



We are attempting to fit partially correlated data including priors for some of 
the parameters. With our method we closely follow the approach discussed in 
ref. [70]. 

The data for / PS , mps, r /a and Zp 1 ' at the various values of the lattice spacing 
and the different volumes is described by the following set of fit parameters 

r */o, r$B , rffA^, k M ,F, C 1)2 , K/a}^ {zf E (2 GeV)} p , D mps , D fps . 

With the notation {...}/3 we mean that there is one fit parameter for each 
lattice spacing involved in the fit. Those parameters can be used, employing 
eqs. (8, 12) and (13), to compute predictions for amps, a/ps, r /a and Z^. 
We shall denote these predictions by 

and they all depend on a sub-list of fit parameters and on the spatial extend 
L/a, as presented in section 4.1. 

We adopt the following definition for the x 2 -function at fixed value of (3: 



(B.l) 



Clflq 



S[r /a] J V S[Z, 



The matrix C is the properly normalised covariance matrix, see e.g. refs. [71, 
72]. The last two terms in eq. (B.l) are included in order to appropriately 
include the measurements for r /a and Zp 1 ' . Values for the correlation among 
mps and fps are given in the tables C.6-C.9. The full x 2 is obtained by sum- 
ming x 2 (P) over & 11 /^-values. 

In case we need to add prior knowledge for a fit parameter 9 in order to 
stabilise the fit we include additional terms of the form 



'Vo-9\ 2 



5[9\ 



(B.2) 



to the fit, assuming Gaussian error distribution for 9 with standard deviation 
S[9}. Here the predicted values Vq are the fit parameters and the "data" with 
errors are the priors. The term eq. (B.2) follows from Bayesian statistics, see 
ref. [70] and references therein. 
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In detail the computation of x 2 is performed in the following steps for fixed 
values of the fit parameters: 

(1) compute infinite volume and continuum predictions: 

compute the continuum predictions for mpg and /pg using the formu- 
lae (8) in infinite volume for given values of rg/iR. Determine r^/a for 
each /3-value according to eq. (12) and use it to estimate the scales 1/a 
from f v . 

(2) compute finite volume predictions in the continuum: 

now apply finite size corrections to scale the data to the finite volumes 
used in the lattice simulations. 

(3) compute finite a predictions: 

use rj/a to compute the estimates of ampg and a/pg for finite values of 
the lattice spacing a, depending on the coefficients D mps j PS . 

(4) compute x 2 '- 

the finite volume and finite a predictions for mpg and /pg can now be 
compared to the data and the x 2 can be computed as detailed above. The 
X 2 includes also the corresponding terms for the chiral extrapolation of 
r /a and for the estimate of Z™ s . 

The errors are estimated using the bootstrap method. For mpg and /pg we use 
bootstrap samples as generated from the raw data. For r /a and Zp l those 
are generated from a Gaussian distribution with mean and standard-deviation 
properly chosen. This is justified since we checked that the correlation of m PS 
and /pg to r /a and Zp 1 is negligible. This procedure allows to estimate errors 
on primary quantities, such as fit parameters for instance, as well as on derived 
quantities, like e.g. the lattice spacing. 
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C Data Tables 





amps 


a/ps 


ampcAC 


rn/a 


LI a 

1 


A x 0.0060 


0.1852(9) 


0.0770(8) 


+0.0019(4) 


4.321(32) 


24 


A 2 0.0080 


0.2085(8) 


0.0835(4) 


+0.0008(3) 


4.440(34) 


24 


A 3 0.0110 


0.2424(5) 


0.0892(3) 


-0.0002(5) 


4.362(21) 


24 


A 4 0.0165 


0.2957(5) 


0.0969(2) 


-0.0017(2) 


4.264(14) 


24 


A 5 0.0060 


0.1831(6) 


0.0784(4) 


+0.0005(4) 


NA 


20 


Table C.l 












Data at (3 = 3.8 












aflq 


tin (<py 


7 

a /ps 


"''tPCAC 


-, 

r o/ a 


T In 

J-J 1 Li 


Bi 0.0040 


0.1362(7) 


0.0646(4) 


+0.00017(17) 


5.196(28) 


24 








— n nnnnQi'1 7^ 




94 


B 3 0.0085 


0.1940 5) 


0.0742(2) 


-0.00052(17) 


5.130(28) 


24 


B 4 0.0100 


0.2100(5) 


0.0759(4) 


-0.00097(26) 


5.143(25) 


24 


5 5 0.0150 


0.2586(7) 


0.0830(3) 


-0.00145(42) 


5.039(24) 


24 


B 6 0.0040 


0.1338(2) 


0.0663(2) 


+0.00022(11) 


5.259(21) 


32 


5 7 0.0030 


0.1167(4) 


0.0633(3) 


+0.00030(14) 


NA 


32 


Table C.2 












Data at /3 = 3.9 












afi q 


ampg 


a/PS 


ampcAC 




r /a 


L/a 


d 0.003 


0.1038(6) 


U . \JO\JKJ I t; 1 


-4-0 nnnsfini") 


6 584f34l 


39 


C 2 0.006 


0.1432(6) 


0.0569 2) 


-0.00004(14) 


6.509(38) 


32 


C 3 0.008 


0.1651(5) 


0.0595(2) 


-0.00065(13) 


6.494(36) 


32 


C 4 0.012 


0.2025(6) 


0.0644(2) 


-0.00092(14) 


6.284(22) 


32 


C 5 0.006 


0.1448(11) 


0.0558(5) 


-0.00027(19) 


NA 


24 


C 6 0.006 


0.1520(15) 


0.0508(5) 


+0.00002(20) 


NA 


20 


Table C.3 












Data at = 4.05 










CL[l q 


ampg 


o/ps 


ampcAC 


r /a 


L/a 


Di 0.0020 


0.0740(3) 


0.0398(2) 


+0.00006(6) 


8.295(45) 


48 


L> 2 0.0065 


0.1326(5) 


0.0465(3) 


-0.00032(11) 


8.008(29) 


32 


Table C.4 












Data at /3 = 4.20 
















' = 1/a) 










3.80 0.340(10) 







3.90 0.384(08) 
4.05 0.417(11) 

Table C.5 

Data for Zp values in the RI'-MOM scheme at scale fi' = l/a. 
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Ai A 2 


^3 A 4 




-0.72 -0.66 


-0.03 0.06 




Table C.6 

Correlation among a/ps and amps at (3 = 


3.8 




B\ B2 B% B L 


1 B 5 Bq 


B 7 


-0.53 -0.72 -0.44 0.10 0.34 -0.62 


-0.74 


Table C.7 

Correlation among a/ps and amps at @ = 


3.9 




Ci c 2 c 3 


c 4 c 5 


Cq 


-0.39 -0.44 -0.34 


-0.07 -0.41 


-0.44 


Table C.8 

Correlation among a/ps and amps at (5 = 


4.05 




Di 


D 2 




-0.44 


-0.22 




Table C.9 

Correlation among a/ps and amps at f5 = 


4.2 




arrip S am ps 


amy 




B x 0.1362(7) 0.109(07) 
B 2 0.1694(4) 0.134(10) 
B 3 0.1940(5) 0.169(11) 
B 6 0.1338(2) 0.110(08) 


0.404(22) 
0.422(09) 
0.428(08) 
0.416(14) 


0.391(15) 
0.434(19) 
0.424(14) 
0.409(21) 


Ci 0.1038(6) 0.090(06) 
C 2 0.1432(6) 0.123(06) 


0.335(12) 
0.347(08) 


0.352(23) 
0.344(13) 


Table CIO 



Comparison of am ps to am ps and m v to amy for the ensembles where the corre- 
sponding neutral masses have been computed. 
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